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We study the vortex lattices in a Bose-Einstein Condensate in a rotating anisotropic harmonic 
trap. We first investigate the single particle wavefunctions obtained by the exact solution of the 
problem and give simple expressions for these wavefunctions in the small anisotropy limit. Depending 
on the strength of the interactions, a few or a large number of vortices can be formed. In the limit 
of many vortices, we calculate the density profile of the cloud and show that the vortex lattice 
stays triangular. We also find that the vortex lattice planes align themselves with the weak axis of 
the external potential. For a small number of vortices, we numerically solve the Gross-Pitaevskii 
equation and find vortex configurations that are very different from the vortex configurations in an 
axisymmetric rotating trap. 
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The study of rotating Bosc-Einstcin Condensates (BEC) have progressed rapidly over the last three years. After 
the initial demonstration of vortices at JILA and ENS regular lattices containing hundreds of vortices have been 

, formed These vortices were found to form remarkably regular triangular lattices. More recently, the dynamical 

properties of BEC's containing vortex lattices have been investigated 0,0]. 
t-H , The success of experimental groups in creating and probing vortex structures have led to a flurry of theoretical 
■ activity. The static and dynamical properties of a BEC with a vortex lattice 0, S IS HH 10 > as wen as vortex lattices 
of spinor and multi component condensates p^, IT^ . IT3 | are active areas of theoretical research. Also, the presence of 
multiply quantized vortices 0, ^| and melting of the vortex lattice |0, are being discussed. The aim of theory 
is to clearly understand the behavior of a BEC with a vortex lattice under different conditions. The versatility of the 
dilute gas BEC experiments enables theorists to consider different, experimentally realizable, scenarios to discuss the 
response of the vortex lattice to various external stimuli. 

In this paper, we consider a BEC in a rotating anisotropic trap. This system is different from most of the exper- 
imental realizations of rotating BEC's to this date . In most vortex lattice experiments, and the theoretical works 
that try to explain them, the rotating gas is first given angular momentum and then placed in an axisymmetric 
trap, which conserves this angular momentum (for an exception see Ref. 19]). The gas then rotates with an angular 
frequency that will minimize its free energy. In this paper, however, we will assume that the BEC is placed in an 
axially anisotropic trap that is rotating with a fixed angular frequency fi, and is in its equilibrium state in this trap. 
Theoretically a similar scenario was considered in Ref. [20l l2l| . however both of these works considered vortex states 
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with small number of vortices and did not present any calculations pertaining to a large vortex lattice. 

Recently such an anisotropic rotating potential was applied to a BEC with a vortex lattice to induce collective 
oscillations 5]. In that experiment, the gas is rotating at a much higher angular frequency than the applied anisotropic 
potential, thus the obtained results are not directly relevant to the theory presented here. It is, however, important 
that a rotating anisotropic trap was demonstrated experimentally and we believe that the system we consider in this 
paper can be realized in a similar experiment. An important result of that experiment was that a vortex lattice that 
is not triangular was observed, although it was not an equilibrium state. One is naturally led to ask whether such 
a structural change in the lattice is possible as an equilibrium state, when the anisotropic potential is static in the 
rotating frame. The main result of this paper is that, when the cloud contains a large number of vortices, the deviation 
from the triangular lattice is small and no structural phase transition should be expected. However, the configuration 
of a small number of vortices can be markedly different from the vortex configurations of an axisymmetric trap, in 
agreement with Ref.j^]. Also, an anisotropic trap was used to study vortex generation [l^, and to investigate the 
escape of a BEC from the trap |22( . The scenario considered in this paper may be experimentally realized in a setup 
similar to one of these three experiments. 

The paper is organized as follows. In section II, we find the single particle wavefunctions for a rotating anisotropic 
trap by exactly diagonalizing the Hamiltonian. We then show that the wavefunctions in the lowest Landau level (LLL) 
can be written as simple analytic functions in the small anisotropy limit. Sections III, IV and V contain discussions 
related to a condensate with a large number of vortices. In section III, we give a variational wavefunction for BEC 
with a vortex lattice in the lowest Landau level of the anisotropic trap. We then use this wavefunction to calculate the 
density profile of the cloud as a function of the rotation frequency and the anisotropy of the trap. In section IV, we 
show that the vortex lattice stays triangular and calculate the small distortion of the lattice caused by the anisotropy 
of the density profile. In section V, we show that the minimum energy configuration of the lattice corresponds to 
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lattice planes aligning with the weak axis of the external potential. Section VI contains a discussion of the vortex 
structures when only a small number of vortices are present. We present the results of the numerical solution of the 
Gross-Pitaevskii equation using the wavefunctions found in section II. Finally, we give a summary of the results in 
section VII and discuss their consequences for experiments. 



II. SINGLE PARTICLE EIGENSTATES 



To be able to discuss the properties of a BEC in a rotating anisotropic trap, we first consider a single particle in 
such a trap. This problem is exactly solvable (see [2l| and references therein) and this solution guides us when we 
introduce a variational wavefunction for the vortex lattice in the next section. Also, we use a truncated basis of the 
exact single particle wavefunctions to study the BEC containing a small number of vortices in section VI. 

We start by introducing the Hamiltonian for a two dimensional axisymmetric trap 



An anisotropic trap will have an additional term 



1 



(1) 



(2) 



In a reference frame rotating with an angular frequency f2', the system is described by adding a term —tt'L to the 
Hamiltonian. Using lengths measured in units of oscillator length 
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and measuring energy in units of Tiujq, the Hamiltonian is 

U=\ [(l + V )x 2 + (l- V )y 2 +pl 



tt{xp y -yp x ). 



Here the dimensionless rotation frequency is 



and the dimensionless anisotropy is 



A 2 



(3) 
(4) 
(5) 
(6) 



In the limit of r\ — > 0, we recover the Hamiltonian for a rotating axisymmetric trap. 
This Hamiltonian can be diagonalized by the linear canonical transformation 

u' = j x (cos(0)x - sia((f))p y ) 
q'u = — (sin(^)y + cos(0)p x ) 

lx 

v = j y (cos((p)y - sm(4>)p x ) 



(7) 



where 



and 



— (sm((j>)x + cos((j))py) , 
ly 



tan(20) = — , 
V 



lx = 



ly = 




(8) 



(9) 
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The diagonal Hamiltonian is then 



H = -Wi + 2 A /r»2 



W (u" + q> u 2 ) 



(10) 



, q' u and are canonically conjugate pairs, the energy spectrum in units of Tiluq is simply 



e = ne„ + me v + -(e u + e„), 



(11) 
where 

e u - \j\ + 2sj& + (|j I +« 2 (12) 

e„ - ^l-2yJw+(ff + W, (13) 

and n, m are non-negative integers. To get a spectrum that is bounded from below, we must have 

T^l-ft 2 . (14) 

If this condition is not satisfied, the centrifugal force overcomes the trapping potential in the weak direction and the 
particle will no longer be trapped. However, when this condition is satisfied as an equality, e v is zero and there is an 
infinite set of degenerate eigenstates forming the analogue of the lowest Landau level (LLL) 23] for an anisotropic 
trap. The formation of such a Landau level will simplify the discussion of the vortex lattice considerably. 

To obtain a useful form for the eigenstates of the anisotropic trap in real space, we express the creation-annihilatk 
operators of the anistropic rotating trap in terms of the creation-annihilation operators of the the axisymmcti 
rotating trap. For the axisymmetric trap we define 
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-j=(a x + ia v ) 



Similarly, for the anisotropic trap, we have 



a u = A=(u' + iq' u ) 



aov = —7={ a y +ia x )- 
V2 



We can express the latter in terms of the former as 



a„ = - 



a„ = — 



{lx + — )oi + {lx - —)a\ 

lx lx 



(15) 



(16) 



(17) 



with 



<2l = 



0-2 = 



1 

V2 
1 

V2 



[(c + s)aou + i(s - c)a 0v ] , 



(18) 



[(c + s)a 0v + i(s — c) 



The ground state of the axisymmetric Hamiltonian, |00)o, satisfies 

a Ou |00) o = a 0v \00) = 0. 

This state is also annihilated by ai,a2 

ai|00) = a 2 |00) = 0. 



(19) 
(20) 
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Similarly, the ground state of the anisotropic rotating trap, 1 00) is defined by 

a u |00) = a t ,|00> = 0. 
Rewriting the definition of a u <|17ll . we have 

1 



with 



{-fx + — )ai + {"f x - —)a\ 

lx lx 



-e 2 



(21) 



(22) 



Qx = 



e„ = 



1 X 

(7. + —) 

lx 

Which gives us an easy way to identify the ground state of the anisotropic potential 

= a u \0) = 2e-^^ 2 ai e^ a ^\Q) . 

By reexpressing both a v and a u as above we obtain the correctly normalized ground state as 



(23) 



(24) 



|00) 



1 W 1 N 

(7, + — ) (ly + — ) 

lx ly 



-\®M? -5®«(4) a lnnN 

e I e £ |00)o- 



(25) 



Higher states in the anisotropic rotating trap are 



\nm) = — = =100), 



and their overlap with the ground state of the axisymmetric trap can be calculated as 



o (00|nm) 



1 



mm 



>n— lnm- 1 



\ 
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1 



(lx + —){ly + —) 



lx 



1 

ly ' 



(26) 



(27) 



These overlaps will help us represent the wavefunctions of the anisotropic potential in a simple form. 

For the axisymmetric rotating trap, wavefunctions in the LLL are analytic functions of z = x + iy. When the 
rotation frequency is close to the trapping frequency, these wavefunctions become the only states that are populated 
due to the high energy cost of the higher Landau levels. These wavefunctions have successfully been used to describe 
the properties of the vortex lattice, even if the experimental situation does not correspond to the fast rotation limit 
0, la Il2j| . The analytical properties of a many-body wavefunction formed in the LLL describes vortex formation very 
well, thus such a wavefunction can always be used as a variational wavefunction at any rotation frequency. Also, very 
recently, a BEC that is in the LLL regime was experimentally realized [24j. 

For an anisotropic rotating trap, if the rotation frequency is close to the critical frequency defined by Ea.Q14[l. 
Landau levels are again formed and we can consider only the wavefunctions in the LLL to describe vortices. However, 
there is no simple expression for the wavefunctions in the LLL of the anisotropic trap except in the limit, 



V_ 

26 



S = 1 -0< 1 

■q < 1 



(28) 
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This limit corresponds to small anisotropy of the trap and fast rotation, and is as much experimentally accessible as 
the fast rotation limit of the axisymmetric trap |tI l24j. It is also important to realize that, even if the anisotropy 
of the trap is small, as the rotation frequency approaches the critical frequency, the cloud of atoms will be stretched 
along the weak axis. This enables us to use the wavefunctions found in the small anisotropy limit to describe a BEC 
with any aspect ratio in the plane of rotation. For a trap with large anisotropy such wavefunctions should only be 
regarded as variational wavefunctions, while they should describe the system very well in the small anisotropy limit. 

Thus, we concentrate on the small anisotropy limit and find an expression for the real space representation of the 
wavefunctions in the LLL. In this limit various parameters in the exact solution reduce to 



cosO) = sin(0) = -j= 



lx = 1 

ly = 



1 - - 

_ 2_ 

1 + A 

x ~ 2 



(29) 



(W) f -(i+ir 



A\2 



Using these parameters, we calculate the spatial wavefunctions. Recalling that the ground state of the axisymmetric 
trap is 



(zy|00) = -^=e 

V7T 



(30) 



with 



The operators take the form 



z = {x + iy)/V2. 

a 0u = (z + d s )/V2 
a 0v = i(z + d z )/V2. 



(31) 
(32) 



By using (fT?T)f> we have 



(xy\00) 



ly + ly 1 



(33) 



We can also find the form of the excited states as in this limit 



il = K t , + y a 0v ). 



(34) 



We have then 



(xy\0n) = 



1 



V ly + ly 



[z-6 y z~(d, + 9 y d z )] n e e y z - zS . 



(35) 



We define 



and use the transformation 



/„ = [z - 9 y z - {d- z + y 8 z )] n e y z - z \ 



w = z — u,,z 



y . 



(36) 
(37) 



to calculate /„ as 



In = (l-e 2 y )^d y 2 H n ( 



ly Uy 



z)e°^ 



(38) 
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Where H n are the Hermite polynomials. 

Thus the wavefunctions in the lowest Landau level are given by 



1 (-»)' 




(39) 



We check the orthogonality of these wavefunctions using a method which demonstrates how matrix elements of 
operators can be calculated in this basis. The overlap of two wavefunctions |\& n ) and |\& m ) is proportional to 



(*n|* ro ) oc / dxdye 



-2zz+9 y (z 2 +z 2 ) 



H n (az)H m (az), 



with 



We can use the generating functions for Hermite polynomials 

d 



H n (az) 



t +2taz 



dt 



and then evaluate the gaussian integral. Which gives us 



(*«|*r. 



d n d rn 
dt n ds r 



t=o 



t=s=0 



(40) 



(41) 



(42) 



(43) 



This expression is equal to zero unless n = m. With the correctly normalized wavefunctions, we obtain the orthonor- 
mality relation 



(44) 



In the small anisotropy limit presented here, we see that the anisotropic single particle wavefunctions are still 
analytic functions of z, multiplied by a Gaussian. Thus, they only have overlaps with the LLL of the axisymmetric 
rotating trap. The introduction of a small anisotropy reorganizes the LLL within itself, i.e. no wavefunctions from 
higher Landau levels of the axisymmetric trap are mixed. 

The n th single particle wavefunction of the axisymmetric trap corresponds to a state that has a single multiply 
quantized vortex z n at the origin. For the anisotropic trap, however, the n th single particle wavefunction will have n 
singly quantized vortices that are distributed on the weak axis of the potential. Positions of these vortices are given 
by the zeros of the Hermite polynomial H n (See Fig[l]). 

The expression for the single particle wavefunctions in the small anisotropy limit Ea. l|39|l . is one of the central 
results of this paper. In the following sections, we will use these wavefunctions to discuss the properties of a vortex 
lattice in an anisotropic rotating trap. 



III. VARIATIONAL WAVEFUNCTION FOR THE VORTEX LATTICE 



After solving the one particle problem exactly, we now consider a BEC containing many particles in a rotating 
anisotropic trap. Here the system is very well described by the Gross-Pitaevskii energy functional 



E = dx dyV* 



~ [(1 + r))x 2 + (1 - V)y 2 +Pl+Pl] - VL Z 



* + | j dx dy\$>\ 



(45) 



In the fast rotation limit considered in this paper, the behavior of the cloud in the z direction can be easily described 
by Thomas-Fermi approximation. The validity of this approximation has been experimentally checked by viewing 
the vortex lattice from the side 0, l20| . Experimental results in an axisymmetric trap confirm that the vortices are 
only marginally bent along the z direction, as predicted by the Thomas-Fermi approximation We also expect 
the vortex lattice in an anisotropic trap to show the same property and have straight vortices through the cloud. 
Therefore in the following discussion we will only consider the Gross-Pitaevskii energy functional in two dimensions. 
This will allow us to use the results of the previous section more directly. 
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FIG. 1: Density of the 5 eigenstate at various anisotropies. The parameter A was chosen as 0.05, 0.35, .65, .95 starting from 
top left figure. The zeros of the density are placed on the weak trapping axis, at the roots of the 5 th Hermite polynomial, as 
can be seen from the form of the wavefunction Ea.(39| . Length is measured in units of oscillator length. 



In the small anisotropy limit introduced in the previous section, we have found that the LLL of the anisotropic trap 
only has an overlap with the LLL of the axisymmetric trap. Thus any wavefunction in the lowest Landau level of the 
anisotropic problem can be written as 

*(as, y) = /(z)e fl »*V" = g(z) e - zg , (46) 

here both f(z) and g(z) are analytic functions of z. The zeros of the function f(z) (or g(z)) correspond to the positions 
of the vortices in the BEC. Thus, a BEC containing a vortex lattice will correspond to an analytic function f(z) that 
has a regular array of zeros. An analytic function that both has a regula r array of zeros and that gives a normalizable 
wavefunction can be introduced using the Jackobi theta function [12, l25l l26| . We write down the following variational 
wavefunction, which we will refer to as the variational vortex lattice wavefunction 

®(x,y) = ce(Ce llfi ,T)e-^ {xsinv+vcosv} (47) 

xe 2^( x sin cos ip) e 2l±L( x 2 -y 2 +2ixy) e -2-^p£ _ 

Here C is a normalization constant, if determines the angle of the lattice basis vectors with the x axis. We define the 
lattice basis vectors a\ and a,2 as complex numbers in the z = x + iy plane and vq is the volume of the unit cell 

v = Im{a\a 2 }. (48) 

Then ( and r are defined as 

C = */M, (49) 

and 



r = u + iv = a-ijax- 



(50) 
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Finally 7 is a variational parameter that is in general complex, and is used to ensure that the wavefunction is 
normalizable. 

The Jackobi Theta function 

0(C,r) = -i £ ( _ 1)V - (n+ i/2)V™c(™+i/2) (51) 
n=— 00 

is a quasi-periodic function. The first line of the variational wavefunction Ea.|47j. the Theta function with the 
exponential factor, is a periodic function in the complex plane p7| . Similar wavefunctions have been used to discuss 
vortex lattices in superconductors [2{| , superfluid He 4 |2y and BEC in axisymmetric traps 0, H, . 

The energy functional takes a simple form in the LLL as the projection of the angular momentum operator L on 
to the LLL is 



dx dy^*{z)L^{z) = J dx dy(r 2 - l)|^(z)| 2 . (52) 

Thus, the Gross-Pitaevskii energy functional takes the form 

E = 5 J dx dy\^\ 2 (x 2 + y 2 ) + | / dx dy\^\ 2 (x 2 - y 2 ) + | / dx dy\^\ 4 . (53) 

We minimize this energy functional with respect to the variational parameters in the wavefunction to understand 
the behavior of a BEC containing a vortex lattice. The first property we arc interested in is the overall density profile 
of the cloud. 

Unlike a BEC in an axisymmetric trap which preserves its rotational symmetry and spreads out evenly in all 
directions, a BEC in an anisotropic trap will prefer to stretch along the weak axis of the trap. Even if the anisotropy 
of the potential is small, if the gas is rotating close to the trapping frequencies the cloud will be stretched along the 
weak axis of the potential. 

To find the effect of the anisotropy on the cloud profile, we can average out the vortices to write an 'averaged vortex' 
wavefunction. This will correspond to replacing the periodic part of our variational wavefunction with a constant 



(54) 



In this approximation, the orientation of the vortex lattice with respect to the potential is not important, as one can 
always choose the variational parameter 7 to compensate for the effect of the orientation of the lattice on the averaged 
density profile. Furthermore, it is evident from the energy functional that cloud will be stretched along the weak axis 
of the potential. We see that the averaged cloud density will be a gaussian of the form. 

* = Ce-i a * x2 -i a * y2 . (55) 
Doing the variational calculation, we find the half width of the cloud in both directions as 



(1 + A) V 2nS 



-1 _ 1 (l-A 2 ) 1 ^ g N 



a y = 



(1 - A) V 2tt5 



and the vortex density Uq as 



7T , 1 2nd 



(57) 



(58) 



2^o (l-A 2 )!/4y gN - 

Here, N is the total number of particles in the cloud and is related to the normalization constant C as 

\C\ 2 = V - -N. (59) 



We see that as rotation frequency approaches critical rotation 5 = t]/2, the cloud is stretched along the y direction 
and becomes more one dimensional. The averaged vortex approximation breaks down when the minor axis of the 
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Aspect ratio of the cloud 




3 4 5 6 7 8 

Normalized rotation frequency 2 (1-£2)/r| 



FIG. 2: The aspect ratio of the cloud is defined as the ratio of the radii of the cloud in x direction and y direction. Aspect 
ratio is plotted against A -1 , which is the scaled rotation frequency. Here A -1 starts from the critical value 1, and increases. As 
the rotation frequency decreases the aspect ratio becomes closer to 1 which is the value for an axisymmetric trap. 



cloud becomes as small as the distance between the vortices. The aspect ratio of the cloud ^ is plotted in Fig. [2] 
as a function of rotation frequency. This aspect ratio can be experimentally measured in a set up such as Ref . |22| . 
without the experimental sophistication needed to image individual vortices. It is worth noting that this aspect ratio 
is smaller than the aspect ratio of the one particle ground state in the same trap. 

We note here that the result obtained in the averaged vortex approximation can also be obtained by an appropriate 
generalization of the formulation in Ref. 0, where an averaged velocity field is used as a hydrodynamic variable. After 
calculating the vortex density and the density profile of the cloud, we consider the structure of the vortex lattice in 
the next two sections. 



IV. STRUCTURE OF THE VORTEX LATTICE 



To be able to discuss the structure of the vortex lattice, we need to go beyond the averaged vortex approximation 
of the previous section. Thus, instead of replacing the periodic part of the wavefunction by a constant, we need to 
consider all the Fourier components. As the energy functional e depends only on the absolute value of the wavefunction 
but not on the phase, we find it useful to write 



Here g{r) is periodic with 



where 



|*(r)| 2 = Cg{f)e 



l(r + R) = g{f), 



R = ndi + mai 



(60) 



(61) 



(62) 



are the lattice vectors of the vortex lattice. Also, a\ and ai are the basis vectors of the unit cell and n,m are integers. 
This periodicity allows us to introduce a Fourier series representation for g(r) 



g(r)= 9o 



1+ E 

* 90 



9_K_ -iKr 



(63) 



Where K are the reciprocal lattice vectors 



K = nb± + m&2, 



(64) 
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with n,m integers and the reciprocal lattice basis vectors given by 



b\ = z x a\ 
62 = &2 x z. 



(65) 



To find the equilibrium vortex lattice configuration we need to minimize the energy functional with respect to the 
lattice basis vectors. In doing this minimization, it is important to notice that if there are many vortices in the cloud, 
there will always be a large number of vortices along both the minor and the major axes of the cloud 



a x v a < 1 

OLyVQ < 1. 

Thus, for all reciprocal lattice vectors K except K = 0, we have 

a x < \K\ 
a y « \K\ 



(66) 



(67) 



This gives us the important result that the contribution of the higher fourier components of the vortex lattice to 
the first two terms in the energy functional are suppressed by the exponential factors e**' 1 ** or e »' » . So, if the 
conditions Eq. (|66|l are satisfied, i.e. the cloud is not extremely stretched, it is sufficient to consider the minimization 
of the interaction term in the energy functional 



dx dy\^\ 



C I dx dye 



-2a±x* 



9R c -iK-f 



(68) 



Once again in the multiplication of the fourier components any terms that depend on r will be exponentially small. 
Thus, the optimum lattice is found by the minimization of 



E 



K 



M 2 



(69) 



This, however, is exactly the minimization done for He 4 and BEC in axisymmetric traps |l2i l26j|. Resulting vortex 
lattice is, as in the other cases, triangular. We then find that the vortex lattice of a BEC containing a large number 
of vortices stays triangular even if the confining potential is anisotropic. 

Physically, the robustness of the vortex lattice can be explained as follows. If there are enough vortices in the 
cloud, the stress created by anisotropic confinement is completely screened and the structure of the vortex lattice is 
determined locally by the interaction between the vortices. Unless the trap is not quadratic, or the cloud is extremely 
stretched, the minimum energy configuration of the vortices will correspond to a triangular lattice. 

It is, however, instructive to calculate the corrections to the triangular lattice. The first correction to the lattice 
will come from the change in the cloud's shape, not directly from the anisotropy of the trap. Thus, the deviation from 
the triangular lattice will increase closer the rotation frequency gets to the weak trapping potential. To calculate this 
deviation, we once again turn our attention to the interaction term in the energy functional. 



ej = ^-C z I dx dye 



-2alx 2 -2a 2 y 2 



V 

K 



9k -ilt-r 



(70) 



Now instead of considering terms that have no position dependence, we need to consider all the terms that have 
position dependence only in the strong confining direction, i.e. the minor axis of the cloud. This is due to the fact 
that as the cloud gets more and more stretched we get 



Thus, to study the structure of the lattice we need to minimize 



(K x —K ' ) 2 



(71) 



(72) 
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Change in the lattice parameters as the cloud is stretched 



0.87 - 
0.8660 - 




FIG. 3: The separation between the lattice planes u (see Ea. l^Ojl fis plotted as a function of o> x vq. This parameter a x uo can 
be regarded as one over the number of vortices along the minor axis of the cloud, u — v/(3)/2 corresponds to a triangular 



lattice, and we see that the deviation from the triangular lattice remains less than 0.5 % even around a x vo 
theory presented here should break down. 



1, where the 



Now, it is important to consider the orientation of the vortex lattice. In the next section we show that the lattice 
planes prefer to orient themselves along the weak trapping direction, thus we choose the reciprocal lattice basis vectors 
as 



&i = b\x 

b 2 = bi(ux + vy). 

The fourier components g K can be calculated by expanding the 6 function 0,H3 as a sum of Gaussians 

9R = M) mi ~ 



(73) 



for 



Thus, we minimize 



K = mi&l + 77J2&2- 



J _ ^(mi + l)(ru-mi)/2 e -^ [(™i) 2 + (™i) 2 +2(in2-nmi) 2 -^(»(m 2 -n2)) 2 ] 

mi ,m2 ,ni 



(74) 



(75) 



(76) 



This minimization can be done numerically, as the exponential dependence causes a very rapid convergence of the 
sum. We find that the triangular lattice becomes compressed in the weak confinement direction as the rotation 
frequency approaches confining frequency. The lattice planes, all of which are aligned in the weak confining direction, 
move apart while keeping vortex density constant. This increases the density of vortices inside each plane, as can be 
expected. We must remark here that this compression is a very small effect, even at o. x vq ~ 1 where the theory is 
expected to break down, the fractional change in the lattice constant is less then 0.5%. It seems very unlikely that 
any deviation from triangular lattice to be seen in an experiment. As a function of ol x vq this dependence is plotted 
in Figure 3. 

We see that unless the cloud is extremely stretched so that there are very few vortices along the minor axis, the 
deviation from the triangular lattice is small. Thus a BEC in equilibrium with a quadratic rotating potential will 
always have triangular vortex lattice and would not show any structural phase transitions. 
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V. ORIENTATION OF THE VORTEX LATTICE 



Upon the formation of a vortex lattice in an axisymmctric trap, the full rotational symmetry of the system is 
reduced to the point symmetry group of the vortex lattice. For a triangular lattice this group is Cq. The angle of the 
vortex lattice with respect to a fixed direction in the rotating frame is, however, arbitrary. This fact is experimentally 
confirmed by viewing the vortex lattice from the side in JILA experiments Q . 

For an anisotropic trap, there is already a fixed direction in the rotating frame that breaks the full rotational 
symmetry. Energy of the vortex lattice depends on the orientation of the lattice planes with respect to the anisotropy 
of the cloud. Using our trial wavefunction we calculate this energy for an arbitrarily oriented triangular vortex lattice. 

The dependence of the energy on the orientation of the vortex lattice will be due to the only non-rotationally 
invariant term in the energy functional 

e 2 = ^Jdxdy\*\ 2 (x 2 ~y 2 ). (77) 

Again we write the cloud density from our variational wavefunction as 

|*| 2 = Cg(r)e- a > 2 -< v2 , (78) 

with g(r) a periodic function 

g(r)=g(R + f), (79) 

and R the lattice vectors of the vortex lattice 

R = nai + ma 2 . (80) 
For a triangular lattice that makes an angle if with the minor axis of the cloud 

di = a(cos(ip)x + sm(ip)y) (81) 
a 2 = a ( ~ (cos(ip)x + sin(</j)y) + (- sin(<^)x + cos(^)y) J , 



and the vortex density is 



v Q = y—. (82) 



By Fourier transforming the density and evaluating the resulting gaussian integrals we calculate the anisotropy 
term in the energy functional as 



(2 



v r 2 n - V" 



n.m 





--) 


( 









2al 2a 2 y 



(83) 



5£ (*W*) 2 a)" 

x 4Qx g 2ay 

9o 

Here the sum n,m is over all positive integers and 



nbi+mb 2 , (84) 



with 61,62 are the reciprocal lattice vectors 



An 

61 = —;=- (a 2 x z) (85) 

V3a 

- 47T _ 

62 = —=- (zx Si). 
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Now the exponential factors get very small as n,m increase, thus by summing over the closest six K points to the 
origin we obtain the first contribution to the orientational energy of the vortex lattice as 



ai V3z/ 



sin ^_^ )e -7S^ sin ^-f) (86 ) 



6 



+ cos 2 (ip)e ^"'"o + sin 2 (<^ + — )e ^<*yo ' 

6 

We see that the lattice planes prefer orienting themselves with the weak anisotropy direction. The physical reason 
for this alignment is that the maximum amount of angular momentum can be gained when there are more vortices 
at the minimum of the confining potential. The situation is then similar to having the first vortices at the center of 
a BEC when vortices are initially formed. 

This pinning of the orientation of the vortex lattice should be the first noticeable effect of the anisotropy of the 
cloud on the vortex lattice. This effect was first pointed out for small number of vortices in Ref.[2(j. The orientation 
of the vortex lattice planes can easily be checked experimentally by viewing the BEC from the side as done in the 
recent JILA experiments. 

VI. BEC WITH A SMALL NUMBER OF VORTICES 

After discussing a BEC containing a large number of vortices, we turn our attention to the case where there are only 
a few vortices in the condensate. Although, at first it may seem that the fast rotation, i.e. LLL, assumption made 
through out the paper would prevent the the discussion of this limit, this is not the case. The number of vortices in 
the condensate is determined not only by the closeness of the rotation frequency to the trapping frequency, but also 
the strength of the interactions between the particles in the condensate [T1II23. A non-interacting gas, for example, 
would not exhibit any vortex states, no matter how fast it is rotating. 

The relevant energy scales for determining the number of vortices are the rotation frequency S, the strength of the 
interactions gN, where N is the total number of particles in the condensate, and for an anisotropic trap the anisotropy 
parameter 77. As the ratio increases from 0, vortices will start appearing in the cloud. So if the interactions are 
weak enough or the number of particles in the condensate are small, a BEC which has a small number of vortices but 
which is rotating fast enough to be entirely in the LLL is possible. Also any vortex state wavefunctions obtained in 
the LLL should be good variational wavefunctions for slower rotating vortex states as LLL wavefunctions adequately 
represent the analytic properties of all the vortex states. 

In the very first experiments that demonstrated the presence of vortices in a BEC, small numbers of vortices 
were observed [l], y, [lj| • These vortices were found to form very regular configurations. We investigated whether the 
regularity of these arrays would be affected if the confining potential were made anisotropic. By numerical calculation, 
we found that the regular configurations of the vortices are sensitive to the anisotropy of the confining potential. Under 
anisotropy, the vortex configurations first oriented themselves to accomodate the maximum number of vortices on the 
weak trapping axis and then, as the anisotropy is increased, became stretched along the same axis. Our approach is 
similiar to Ref . |2l| , and we obtain similar results with our approximate wavefunctions Eq. (|39|l . This can be regarded 
as a success of the approximation scheme introduced in Section II. 

To determine the vortex states, we numerically solved the Gross-Pitaevskii equation using a truncated set of single 
particle wavefunctions found in Section II. We expand the condensate wavefunction ip(r) as 

n=n cu t 

^(r)= C «*«W- (87) 

n=0 

Using this expansion in the expression for the Gross-Pitaevskii energy functional we obtain 

e = ^ e «l c «| 2+ R n mp qC* n c* m c p c q . (88) 

n n,m : p,q 

Here e n are the eigenvalues of the n th single particle state calculated in Section II, and 

Rnmpq = j dx d V*®*n®m^P®q, ( 89 ) 

are calculated using the generating function for Hermite polynomials as 

2 _( n+m+p+g )/ 2 / -l\-2 

/? _ v ' v Z y \0 |(«+™+p+<2)/2 (qo) 

^nrnpq — , — | 1 — |— T \Vy\ W u / 

iry'nimip'.q] 



14 



-4 -2 2 4 




FIG. 4: Density of various 4 vortex states obtained by numerical calculation. The parameter A was chosen as 0.05, 0.25, .85, .90 
starting from top left figure. Interactions were chosen to result in 4 vortex states at each different anisotropy. Length is 
measured in units of oscillator length. 



x {d r ) n {d s ) m {d t y{d u ye 



We numerically minimized the energy functional by varying c n and plotted the density of some of the representitive 
vortex configurations in Fig. [3]. In this minimization the parameters of the Hamiltonian were chosen so that the cutoff 
imposed on the wavefunctions did not significantly alter the results of the minimization. 

As a simple example to demonstrate the effect of anisotropy consider a BEC with four vortices (See Fig. [4]). If this 
state was created in an axisymmetric trap, the vortices form a square, with the origin as the square's center. If an 
anisotropy is introduced, we find that, first two of the vortices align themselves with the weak trapping direction and 
then move a little bit further away from the origin. The other two approach the weak trapping axis, thus, turning 
the square into a rhombus. As anisotropy is increased, the rhombus becomes more compressed, while conserving its 
form. After a critical value of the anisotropy, the vortices, instead of forming a rhombus, become aligned in the weak 
trapping direction. This critical value, however depends on the ratio of interaction strength to the rotation frequency 

gN_ 
5 ' 

A similar scenario plays out for vortex states with higher numbers of vortices in them. The initial highly regular 
configuration becomes compressed in the strong confinement direction and at some critical values of the anisotropy 
one or more vortices jump to the weak trapping axis. In the situation of extreme anisotropy, all the vortices become 
aligned on the weak trapping axis, forming a one dimensional vortex lattice. It is important to notice that, the 
behavior obtained in numerical calculations is in agreement with the analytical calculation made in the previous 
sections. 

In Fig. [4] , we display some of the obtained vortex states, by plotting the two dimensional particle density. We have 
calculated states with up to 8 vortices in them and believe that these configurations can be experimentally observed 
by imposing a rotating anisotropic potential on the few vortex states obtained in experiments 



perimc 
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VII. CONCLUSION 



We studied the properties of vortex arrays in a BEC that is confined by a rotating anisotropic potential. All our 
calculations were made in the experimentally accessible limit of fast rotation, however we argued that results obtained 
under this assumption should also apply to slower rotating condensates. 

To facilitate the discussion of a condensate in an anisotropic rotating potential we first solved the one particle 
problem exactly and obtained analytical expressions for the single particle wavefunctions. 

We obtained a variational wavcfunction for a BEC containing a vortex lattice using the results of the single particle 
problem. This variational wavefunction enabled us to calculate the density profile of the cloud, which could be easily 
measured experimentally. Also, we found that even under anisotropic confinement the vortex lattice stays triangular. 
However orientation of the lattice planes is controlled by the external potential. 

For a condensate with a small number of vortices, we numerically calculated the wavefunction using a truncated 
basis of single particle eigenstates. We showed that the regularity of the vortex configurations observed in an axisym- 
metric trap are disturbed by the anisotropy of the confinement. We found that at high anisotropy vortices form one 
dimensional arrays. We believe these new vortex configurations can be experimentally obtained. 
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